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Abstract 

Random first order transition theory is used to determine the role of attractive and repulsive 
interactions in the dynamics of supercooled liquids. Self-consistent phonon theory, an approximate 
mean field treatment consistent with random first order transition theory, is used to treat individual 
glassy configurations, while the liquid phase is treated using common liquid state approximations. 
Free energies are calculated using liquid state perturbation theory. The transition temperature 
T4, the temperature where the onset of activated behavior is predicted by mean field theory, the 
lower crossover temperature T* where barrierless motions actually occur through fractal or stringy 
motions (corresponding to the phenomenological mode coupling transition temperature), and T^-, 
the Kauzmann temperature (corresponding to an extrapolated entropy crisis), are calculated in 
addition to T*, the glass transition temperature that corresponds to laboratory cooling rates. Re- 
lationships between these quantities agree well with existing experimental and simulation data on 
van der Waals liquids. Both the isobaric and isochoric behavior in the supercooled regime are 
studied, providing results for ACy and AC P that can be used to calculate the fragility as a func- 
tion of density and pressure, respectively. The predicted variations in the a-relaxation time with 
temperature and density conform to the empirical density-temperature scaling relations found by 
Casalini and Roland. We thereby demonstrate the microscopic origin of their observations. Finally, 
the relationship first suggested by Sastry between the spinodal temperature and the Kauzmann 
temperatures, as a function of density, is examined. The present microscopic calculations support 
the existence of an intersection of these two temperatures at sufficiently low temperatures. 
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1. INTRODUCTION 



Owing to their mutual repulsion, molecules in a dense supercooled liquid can only rear- 
range by performing an intricately correlated dance. The growing complexity of the coop- 
erative motion as the liquid cools leads to the glass transition. In the deeply supercooled 
state, a sufficiently large number of molecules move in the cooperative activated events so 
that the dynamics can be predicted rather accurately from the mesoscale thermodynamics 
of the configurations. The microscopic theory of glasses based on random first order transi- 
tions quantitatively predicts^ the typical relaxation time^, its distribution^ and the size and 
geometry^ of the independently rearranging regions using only the knowledge of the config- 
urational entropy which can be obtained to reasonable accuracy from experiment! 5 1 6 1 7 1 8 1 . 

Furthermore the random first order transition theory provides a route to calculate the 
relevant mesoscale configurational thermodynamics directly from the intermolecular forces. 
The theory thus allows one to determine what features of the intermolecular forces are impor- 
tant for the configurational thermodynamics and therefore what determines microscopically 
the glassy dynamics. These topics are the focus of this paper. 

Not unexpectedly the details of the repulsive forces are key. For glass physics a purely 
steric picture envisioning infinitely hard repulsion provides a good place to start, but, as 
we shall see, a truly quantitative treatment must go further. Since the time of Rice^ and 
KirkwoocP^ we have become used to the dominance of hardcore sterics in determining gross 
liquid structure. Alder and Wainwright's computer calculations^ and Widom's^ analytical 
treatment of freezing have made Newton's conjecture on the steric origin of crystals^ a 
quantitatively useful paradigm. Likewise the random first order transition theory shows the 
configurational entropy density of a hard sphere system will fall with decreasing volume 
leading to slowing of the dynamics and eventually results in a glass transition. 

The impenetrable hard sphere picture, we shall see, however, is not quite accurate enough 
to describe real molecular liquids in the supercooled regime. High-pressure experiments make 
this failing clear. Apart from a trivial multiplicative temperature scaling of the time scale, 
the dynamics of hard spheres must only depend on density. In the supercritical and ordinary 
temperature regime above the melting point studies of transport show that diffusion and 
rotational diffusion do indeed depend primarily on density for van der Waal's liquids. At 
constant volume the temperature dependence of transport coefficients is almost negligible^!. 
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Such a weak temperature dependence at constant volume is not found in the glassy regime, 
in contrast. While the apparent activation energy for flow at the glass transition is smaller 
at constant volume than along lines of constant pressure, it is not negligible: the ratio of the 
constant volume and constant pressure activation energies is measured^ to be about 1/2. 
Developing an empirical scaling of transport in this regime has attracted attention, most 
notably in a series of works by Casalini and Roland^- 5 * 16 * 17 * 18 * 1 ^ . The shape of the repulsive 
and attractive forces clearly must matter in the supercooled regime. 

We will show how the random first order transition theory clarifies the origin of this 
isochoric temperature dependence of the relaxation rate. First it helps to recall that the 
validity of an essentially purely steric picture of gross liquid structure and of nonactivated 
transport under usual thermodynamic conditions can be explained by the existence of a 
small parameter: the ratio of a steepness length of the potential to the particle size itself. 
Arguments from liquid state perturbation theory show structural quantities on the scale of 
the particle diameter can be obtained by a so-called "blip function" expansion! 20 * 21 !. The blip 
function, measuring the deviation between the model impenetrable particles and real ones, is 
nonzero only near the turning point of molecular collisions and covers the range of a steepness 
length. The steepness length measures the range of distances a particle can penetrate with 
a typical thermal kinetic energy, and is the order of 1/10 the particle diameter. Since the 
steepness length is small on the particle size scale, treating the particles as impenetrable 
becomes a good approximation. Random first order transition theory shows that another 
small parameter involving lengths emerges in the supercooled regime: the Lindemann ratio 
characterizing the extent of the local caged motions. 

In the supercooled regime, molecules spend most of their time vibrating about fiducial 
locations characterizing the aperiodic free energy maxima. These vibrations also cover a 
length about 1/10 of the particle spacing revealed by both neutron scattering experiments^ 
and theoretical calculations^ 2 ^. 

The smallness of the Lindemann ratio ultimately explains the quantitative successes of 
random first order transition theory (which is an expansion about mean field theory) in 
predicting the universal patterns of fragility in structural glasses and supercooled liquids^. 
The ratio between the Lindemann length and the steepness length is clearly, however, not 
particularly small. Because of this fact, the configurational entropy and the configurational 
heat capacity, which enter the supercooled dynamics, depend on volume and temperature 
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jointly. This joint dependence explains the significant size of the isochoric activation energy 
for relaxation and transport. 

To make these ideas clearer and more explicit in this paper we present detailed calcula- 
tions of the configurational thermodynamics of the supercooled Lennard- Jones fluid using 
the tools provided by the random first order transition theory. These calculations are based 
on density functiona l 24 * 26 * 27 * and self-consistent phonon approache d 23 * 28 * which can be justi- 
fied by elegant replica methodolog y ^ 29 * 30 * 31 * 32 * 33 * 34 * 35 * 36 *. Related theories have also recently 
appeared.* 37 * 38 * 39 * 4 -^ In addition our calculations rely on reasonably clear, but somewhat un- 
controlled, liquid state approximations like those commonly employed for liquids that are not 
particularly dense or highly supercooled. It will turn out these approximations reproduce 
quite well the laboratory observations on the pressure dependence of the glass transition of 
molecular liquids. While these approximations predict an ideal glass transition at low tem- 
perature and high density, in keeping with extrapolated data they cannot address the strict 
existence of an entropy crisis with controlled rigor. Fortunately taking an agnostic attitude 
towards that ideal transition does not impair the application of the present scheme to the 
currently accessible experimental regime. It appears to be impossible to reach the long time 
scales where corrections to the mean field configurational entropy should appear 25 . At the 
lowest density studied, the entropy crisis occurs in a temperature-density range that is close 
to the spinodal temperature for a liquid-gas transition. As suggested by Sastry^, there is 
the possibility of a temperature and density that is both the spinodal temperature-density 
and the entropy crisis temperature and density. We investigate this possibility using our 
theory. 

The organization of the paper is as follows. We first describe the liquid state approxi- 
mations we use to describe the overall thermodynamics of the dense liquid states. We then 
discuss how the self-consistent phonon calculations for a single aperiodic free energy min- 
imum can be carried out for smooth repulsive forces. The configurational thermodynamic 
properties follow by subtraction. We then present the resulting predictions for the Kauz- 
mann temperature, laboratory glass transition temperature, mean field dynamical transition 
temperature, and the laboratory crossover transition temperature as functions of density. 
In the following section we compare the predictions starting from the intermolecular forces 
with the patterns of laboratory experiments. We show how the theory reproduces the known 
Casalini- Roland scaling of the dynamics over the measured range but suggests the possibility 
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of observing deviations from their scaling if a broader thermodynamic range can be probed. 
The theory also correctly predicts the ratio of isochoric and isobaric activation energies at 
the glass transition. We summarize our results in the conclusion. 

2. THEORY 

Microscopic Theories of the Glass and Liquid 

The universal features of the dynamics of supercooled liquids and glasses are explained 
by the random first order transition theory^* 3 * 30 * 42 *. At the mean field theory level the un- 
derlying microscopic framework provides a description of several characteristic transition 
temperatures that are expected in glassy systems. This theory predicts there is a temper- 
ature, T^, at which a dynamic transition occurs, corresponding to the onset of activated 
motion in a liquid. At the mean field level, the mode coupling theory temperature Tl ICT at 
which the density-density correlation function no longer decays^ and coincides with T* A . We 
note, however, that phenomenological fits to mode coupling theory consistently give a lower 
temperature than the ab initio mean field predictions. At the mean field level there is also a 
temperature, T£, which corresponds to an entropy crisis predicted long ago by KauzmannPl 
Random first order transition theory goes beyond mean field theories to analyze the creation 
of so-called "entropic droplets" which describe the nucleation of small regions containing a 
multiplicity of states in an aperiodic free energy minimum. These droplets arise as correc- 
tions to the mean field with a large ensemble of local minima of the mean field free energy 
functional for the inside of the droplet in an initial mean field solution. Their multiplicity 
corresponds to the "configurational entropy". The size of the critically rearranging droplet 
is thus determined by a balance of the configurational entropy and the free energy penalty 
required to create a surface between two distinct aperiodic structures. 

Dynamic heterogeneities consistent with the mosaic structures arising from the entropic 
droplets predicted by random first order transition theory nearly twenty years ago have 
recently been observed directly in supercooled liquids^ 4 -* 45 * 46 ! More recently, random first 
order transition calculations have quantitatively predicted the magnitude of the barriers to 
reorganization of a mosaic region, as well as the fluctuations in barrier height from one 
region to another with using adjustable parameters^. The theory has also been applied to 



6 



the diffusion of probe molecules in o-terphenyP^ and to aging in structural glasses®!, again 
in good agreement with observation. Of particular relevance to our work is the microscopic 
prediction of a typical relaxation time in the equilibrium supercooled liquid 

r = r ° exp U^J (1) 

AFt 32 (2) 



k B T S c /Nk B 

32T K 



(AC p /Nk B ) (T - T K ) 
DT K 



T — T; 



K 



(3) 
(4) 



We see the predictions conform to the well-known empirical relation known as the Vogel- 
Tamann-Fucher law for the a-relaxation time r. In this expression, S c is the configurational 
entropy per spherical unit, called a "bead." At constant pressure, its temperature derivative 
is proportional to AC P , the heat capacity jump at constant pressure. D is often called 
the liquid's fragility. Eqn. [T] can be used to relate the the configurational entropy to the 
relaxation time scale in the liquid: 

,S * 32 (5) 



Nk B In r/r 

For the one hour time scale we see the configurational entropy per bead is predicted to 
be a universal value, S c /Nk B = 0.82. The prediction that D = 32/ (AC p /Nk B ) conforms 
very well to the experimental data for a variety of glasses including both molecular and 
ionic glasses. To make this fit quantitative, one uses the melting entropy to count beads, 
as described by Lubchenko and Wolynes®' and by Stevenson and Wolynes^. A survey of 
experimental data reveals that the predicted microscopic coefficient in this expression of 32 
is only in error by roughly 10%. 

A convenient mean field treatment of the localized phase of a random first order transition 
is provided by self-consistent phonon theory, which envisions the atoms in a glass vibrating 
about fiducial sites in an effective potential due to the other atoms. The effective vibrational 
frequency is calculated self-consistently as a function of density. For the hard sphere glass 
there exists a mean field dynamical transition density pa below which the only solution to the 
self-consistent equations is a vibrational frequency of zero, corresponding to a mobile liquid 
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where the dynamics are essentially those of a renormalized gas; above this density a new, 
non-zero, frequency appears. For short times, the system is therefore dynamically a solid and 
corresponds to a glassy local minimum in the free energy. The difference F g i ass —Fu quid = T S c 
is accounted for by the configurational entropy of configurations. At a high enough density 
the entropy is predicted by the mean field theory to no longer to be extensive. The limit 
S c = defines px- The glass density p g can be calculated but depends on the experimental 
time scale of interest (see below). For systems with potential energy functions other than 
purely hard sphere, there are corresponding temperatures Ta, T g , and T K . While the mean 
field transition to rigidity occurs at T^, fractal or string-like excitations still provide a 
downhill escape from minima at that temperature. Only at a crossover temperature T* 
lower than the mean field T\ do such string motions become activated thus explaining the 
overestimation of T* by ab initio mode coupling. Indeed, such a crossover is found where 
phenomenological fits to mode coupling theory show a transition^ and deviations from 
Vogel-Tamman-Fulcher behavior are seen in a Stickel plot 51 . Thus, microscopic calculations 
based on random first order transition theory predictPa crossover temperature T c at a critical 
entropy 

S c (T c )/Nk B = 1.28 (6) 

T c is the temperature at which a complete transition from collisional to activated transport 
will be noticed. As a consequence of the universality of the crossover entropy, the relaxation 
time at the crossover is predicted to be universal, as indeed is observecP. 

Originally developed by Fixman to study periodic crystalline solids^, self-consistent 
phonon theory was used long ago by Stoessel and Wolynes^ to calculate pa in a monatomic 
hard sphere glass. More recently, Hall and Wolynes^ used this approach to study a model 
of a network glass in which the are both steric repulsions and athermal bonding constraints. 
Those calculations predict the well-known decrease in fragility with increasing connectivity 
empirically observed in mixtures of network formers. We see that the mean field random 
first order transition theory based on self-consistent phonon ideas is not restricted to hard 
sphere potential functions and can therefore be applied directly to the Lennard- Jones glass. 

Microscopic treatments of liquids traditionally take the view that because the structure 
of the liquid is dominated by repulsive forces, a good approach is to develop perturbation 
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theories based on a hard sphere reference system. The main issue is to find the optimal 
parameters for the reference system. At densities less than the glass transition density, such 
approaches provide reasonably accurate overall thermodynamics but at the higher densities 
accompanying the glass transition, less has been explored about the accuracy of the different 
perturbation approaches. In light of the similarity between the "blip" and Lindemann 
length scales because we lack a non-perturbative approach to the Lennard- Jones liquid, it 
is necessary to examine two related but distinct detailed perturbation theories to treat the 
liquid state thermodynamics. For consistency, we also employ a perturbative self-consistent 
phonon theory. We now discuss the perturbation theories for liquid and glass. 

Liquid State Thermodynamics and the Potential 

We carried out the calculations using two perturbation theories for the free energy of the 
Lennard- Jones liquid that relate the properties of the Lennard- Jones system to the purely 
steric hard sphere system. The first perturbation approach separates the potential along 
the lines of the well-known Weeks- Chandler- Andersen (WCA) approximation^^!!. This 
approximation is often used to treat molecular liquids much above the melting point. A 
second approximation, known to be more accurate for high density fluids, was developed 
by Kang, Ree, and Ree (KRRj^El j n both schemes, the pair-potential for the system is 
written as a sum of an attractive and a repulsive term: 

V(r) = v(r) + w{r) (7) 
« v HS (r) + w(r) (8) 

where vhs is a hard-sphere potential whose diameter d is determined by 

J dr[exp(-v(r)/k B T) - exp(-v H s/k B T)}yHs(r) = (9) 

With either separation the Helmholtz free energy is obtained by adding a first order pertur- 
bation to the free energy of the corresponding hard sphere system: 

A » A HS +^ I drg{r)w{r) (10) 
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The pair correlation can also be obtained: 



g(r) « y HS (r)exp(-v(r)/k B T) 



(11) 



where yus is the cavity distribution. The pressure can be obtained by numerical differenti- 
ating the free energy: 



a\ 12 /(j\ 6 
r J V r 



P W™ 



The specific separation of the Lennard- Jones potential 

V(r) = 4e 

due to WCA breaks the potential at its minimum 

v(r) 

w(r) 



V(r) + e r<2 1 /V 

r > 2 1 / 6 ( t 

-e r < 2 1 /V 

V(r) r > 2 1 /V 



(12) 



(13) 



(14) 



The KRR separation uses an optimized but density dependent cutoff for separating the 
potential 



v(r) 
w(r) 



V(r) - F(r) 


r < A 





r > A 


F(r) 


r < A 


V(r) 


r > A 



(15) 



F(r) = V(X)-V'{X)(X-r) 



(16) 
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= r m P < Pi 

\( n \ _ r I <L (P2)(p-pi) 3 [fc(4p 2 -pi-3p)+3(p-pi)] ^ ^ /i 7 \ 

A(pj-r m + i2( P2 - Pl ) a P1SPSP2 U'J 

= a(p) P2 < P 

where r m = 2 1 /V, a(p) = 2 1 / 6 /p 1 / 3 , a'(p 2 ) = (<9a(p)/<9p) p=P2 , fc = 3.026462, pi = 0.99p c , 
p 2 = 1.01p c , and p c o" 3 = y/2/r^. By design, the KRR separation reduces to the WCA 
separation at sufficiently low density. For the remainder of this paper, we use units where 
a = 1. 

For liquids above the melting point, gns and yns are well described by the Percus-Yevick 
approximation. For the high densities considered in this work, it is necessary to employ 
a more accurate modification of these pair functions similar to that made by Verlet and 
Weis. The Verlet- Weis correction improves the behavior of g(r) near contact and removes 
oscillations in g(r) for large r. We follow the procedures of Robles and Lopez de HarcP^ and 
Verlet and Weis^. First, for r > d, we write g(r) as 

g HS (r,d,T]) = g PY (r/d',r)') + 8g{r) (18) 
= gpyij /d! ,r]') H exp(— p(r — d)) cos(p(r — d)) (19) 

where r\ = p*d 3 /Q, rj' = rj — r] 2 /16, rj' = p*d' 3 , gpy is the Percus-Yevick hard-sphere g(r), 
and B and p are to be determined. g(r) has structure extending for several a at the low 
temperatures and high densities studied in this work. Therefore, we used the analytical form 
of Smith and Henderson^ to generate values for gpy(r) out to 8 a. 

For r < d, we follow WCA^and extrapolate ln5g(r) quadratically about r = d: 

\n5 9< (r) = In | - Q + p) (r - d) + \ f 1 - p 2 ) (r - df (20) 

leading to 

y H s(r,d,r}) = y PY (r/d',r]') + 5g < (r) (21) 



Eqn. 21 is used in Eqn. Mto determine the hard-sphere diameter d. 
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B is determined by the value of g{r) at contact, through the pressure equation: 



Pp/p = Z = l + 4r)g HS (d,d,r)) (22) 
= l + ^(gPY(d/d',r ] ') + I pj (23) 

where Z is the compressibility factor, which is determined given an equation of state. Thus, 
given an equation of state B can be found. p is determined by requiring consistency with 
the compressibility equation: 

k B T [^j = l + pjdv (g HS (r, d, V ) - 1) (24) 



or 



POO 

X (Z) = 1 + Anp dr r 2 (g HS (r, d, rj) - 1) (26) 
Jd 

pOQ POO 

= 1 + Aixp \ dr r 2 (g PY (r /d' \rf) - 1) + 4tt / dr r 2 5g(r) (27) 
Jd Jd 

pOO rd' 

= 1 + 4vrp / dr r 2 (g PY (r/d', r/) - 1) + Anp / dr r 2 (g PY (r/d', r,') - 1) 
Jd' Jd 

POO 

+4tt / drr 2 5g(r) (28) 

Jd 



d' 



oc 



= Xpy + 47rp / dr r 2 (g PY (r/d' \rj') - 1) + 4tt / dr r 2 5g(r) (29) 
Jd Jd 

with 

_ l+r? + ?7 2 , , 

(1-77) 

Thus, given an expression for Z, x can be found. Therefore one also finds expressions 
for B and p. We have investigated the thermodynamics that follows from three different 
approximations for the compressibility factor of the steric system. These are first, the 
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Carnahan-Starling expression 5 -^ 1 for Z 



Zcs = {1 _ 7])3 , (31) 



second, a high density expression for Z, due to Mulero, Galan, and Cuadros^ 

Zmgc = 7 [s + 1.96192s 2 + 0.55927s 3 - 1.10721s 4 + 0.55626s 5 (32) 
-0.11923s 6 + 0.00954s 7 ] (33) 



^ 1-e (34) 

e = (35) 



and finally the Salsburg-WoocP^ expression 



PRCP/P 



+ 1 (36) 



with the random closed packed density, P*rcpi = l-2165£p^l. While the Carnahan-Starling 
equation of state gives accurate values for the free energy at moderately high densities, this 
is not true for the other two equations of state used in this work. Therefore, we added a 
small correction to the free energy in order to bring our calculated values of the free energy in 
agreement with calculations tabulated by Carnahan and Starling^. Denoting the correction 
5f, we used 5f = for Zcs, 8f = +.073 iV ks for Zmgc-, and 5f = +0.413 iV ks for Z$w- 
The expressions for Z also give the free energy for the liquid via 



^ = f liq = InpA 3 - 1 + jf (Z - 1) ^ + P - J dvg(r)Pw(r) + 5f (37) 



Thermodynamics of a Single Glassy Configuration 

We use self-consistent phonon theory to describe the free energy of an individual glassy 
configuration^ 3 ^. This theory relies on the fact that the time averaged density can be 
written as a sum of gaussians representing vibrations of atoms about fiducial sites: 

= E (f ^^(-^(r - ^) 2 ) (38) 
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In the present work, the force constants «j are all taken to be equal, but this approximation 
need not be made, if snapshots of liquid configurations are available. In the independent 
oscillator version of the theory, the effective interaction between two atoms is given by 



Making a Taylor series expansion for the effective interaction, a is found self-consistently 
using the relation 



where f2 is the volume of a cell containing an atom, a = is always a solution to the 
self-consistent equation corresponding to a uniform liquid. Above a critical density pa, 
the self-consistent equation gives a ^ solutions that correspond to a glassy state. For 



without the Taylor series expansion for an FCC solid. Without making the Taylor series 
expansion, Ree found a single solution with a value of a nearly zero (a "liquid-like" value) 
at low densities and the appearance of a second solution with a value of a characteristic of 
a localized state (a "solid-like" value) above a density of ~ 1.025. At even higher values of 
the density the "liquid-like" value disappeared. Ree interpreted the appearance of the larger 
value of a as the onset of a stable solid phase. With the cell constraint, when the Taylor 
series expansion is made only a single solution is found that continuously increases with 
density from a liquid-like value to a large solid-like value. The transition from liquid-like to 
solid-like value occurred at approximately the same density as the appearance of the solid- 
like value when solving the complete equations. We note that the density 1.025 is close to 
the values often quoted from simulation for the glass transition density (0.99-1.0 1 65 , 1.0£p^). 
In the present work, we found the Taylor series approximation convenient and therefore took 
the value for p* A = pac 3 (= 1-025) without the Taylor series expansion as the value we used 





(40) 



periodic crystalline solids, Ree^ compared the uncoupled oscillator approximation with and 



for p\. 
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The free energy from the self-consistent phonon approximation is given by 

fus 9 ia SS = = p f dRg{R)0V°" {R) + |ln ^ - 3 In [S(v^)] ) (41) 

D = P — (42) 

where Q = D 3 and we have assumed a cubic cell for convenience. Barker demonstrated^ 
that the self-consistent phonon approach produces an error of -0.723 iV ks (-0.224 N ks) in 
the entropy of the hard sphere crystal near FCC close packing using the uncoupled (coupled) 
oscillator approximation. For glasses, we expect the error in the entropy to be less than the 
one found for the FCC lattice and therefore we use the somewhat smaller correction to the 
coupled oscillator value for the glass. This correction has a magnitude of 0.224 N ks- Thus, 
the expression we use for the free energy of an individual Lennard- Jones glassy configuration 
is 

f glass = P J dRg(R){3V e "(R) + \^^~ - 3 In [$(yfeD)]) - 0.224 (43) 
+ ^ J drg(r)(3w(r) (44) 

With this correction, the theory yields values of p* K ranging from 1.16 to 1.24, depending 
on the equation of state used to determine g(r). For comparison, using a liquid structure 
based density functional theory, Singh, Stoessel and Wolynes^ found p* K = 1.14. 



Pressure, Configurational Entropy, and Heat Capacities 

Random first order transition theory identifies the configurational entropy with the differ- 
ence between the free energies of the glass and the liquid. Thus, the configurational entropy 
is 

TS C = Agi ass -^liquid (45) 
SJN k B = f gl ass - f u quid = Af (46) 

With the present liquid state approximation, the best fits of Af versus temperature were 
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obtained using the function 



A/ = a + a 1 lnT + ^ + || + £ (47) 



AC P and ACy were then be found by differentiation: 



AC.V/.YA7, T[ 9 4l) (18) 



<9T 



AC p /iV ks = T (^j (49) 

The fitting function has the flexibility to describe heat capacities that depend on temperature 
as either 1/T or 1/T 2 ; both functional dependencies have been suggeste d 68 * 69 ^ in empirical and 
model studies. Within the liquid state theoretic incarnation of random first order transition 
theory, we see that no physical significance can be ascribed to these fits such as " Gaussianity 
of the landscape" or inferences of "excitation" structure as in some picturesque models. The 
temperature dependence simply reflects the softening of the steric potentials with increased 
kinetic energy of the molecules. 

The pressure can be evaluated by numerically differentiating the liquid free energy: 

P = P 2 (if) k * T ( 5 °) 



dp 



3. RESULTS FROM SIMULATIONS 



To orient the reader to the present study, we first describe computer simulations of sim- 
ple glasses which have previously been carried out. Computer simulations are often used to 
discuss the structure and dynamics of glassy materials. Yet, despite their widespread use, 
carrying out meaningful simulations of these systems is not easy. Direct studies are hindered 
by several factors. One difficulty is the need to avoid crystallization. The other challenge is 
to reach the very long time scales consistent with laboratory experiments (microsecond to 
seconds). So far this has never been done. For monatomic glasses, even on a short time scale 
avoiding crystallization is extremely difficult. This is one of the reasons binary Lennard- 
Jones mixtures have become popular choices for simulation. When crystallization is avoided 
it is still important to recognize that T g depends on cooling rate, albeit logarithmically. 
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Jonsson and Andersen^ found T* = ksTg/t = .38 in a constant pressure simulation of an 
80/20 mixture (with parameters chosen similar to those used by Stillinger and Weber to 
study amorphous Ni 80 P 2 ^) with a cooling rate of AT* / At* = 0.001. Kob and Andersen^ 
used smaller cooling rates (as slow as 1.5 x 10~ 7 ), fit the diffusion constant to the predic- 
tions of mode coupling theory, and determined T^ CT = 0.435. Sastry, Debenedetti, and 
Stillinger^ used cooling rates ranging from 10~ 4 to 10~ 7 and found the onset of activated 
dynamics at T* < 0.45 at a density p* = 1.2. Sastry fit the temperature dependence of the 
diffusion constant calculated in simulations. The fits give T£ = 0.30 at p* = 1.2, similar to 
values found by Sciorino, Kob, and Tartaglia from an analysis of inherent structures^ and 
by Coluzzi, Parisi, and VerrocchicJ^, who used a replica-based approach. 

Many simulations study monatomic Lennard- Jones systems but modify the Lennard- 
Jones potential by adding a small many-body potential energy term that discourages 
crystallization^. The modified potential has triple point and critical point values T 4 * p = 
0.687, p* tp = Pt. P .(r 3 /e = 0.0030619, p* p = 0.67, T* = 1.16, p* = 0.109, and p* = 0.247 251 , 
which can be compared to the values for the full LJ potential, T t * p = 0.67, p% t = 0.86, T* = 
1.25 - 1.58, p* = .109 - .303, and p* = 0.26 - 0.4CP. Fits of the simulated diffusion con- 
stant to the results of mode coupling phenomenology, (T* — T^ fCT ) 7 give Tlj CT = 0.475 at 
p* = l.G 7 A Calculations using the modified Lennard- Jones potential and a cooling rate of 
4.2 x 10" 4 found 0.3 < T* < 1.13 for densities 0.85 < p* < 1.25^. 

Clearly the cooling rates for the present day simulations are still much faster than those 
used in experiments. According to Angel™, laboratory cooling rates of ~ 0.17 K/s give 
the experimental glass transition at the temperature at which the structural relaxation 
time equals 100 seconds. With some cleverness cooling rates as fast as 10 6 K/s have been 
reached for melt-spun glasses^. In dimensionless simulation units cooling rates as small as 
10- 7 = dT*/dt* and as large as 10" 3 = dT*/dt* have been used in Lennard- Jones simulations, 

1 /2 

where t* = t (-^) • Rescaling the units for argon, using e/ks = 120 K and a = 3.4 A, 
the range of simulated cooling rates is then 10 7 < AT* /At* < 10 11 K/s. For the monatomic 
Lennard- Jones simulation^, therefore, the cooling rate is roughly 10 4 times faster than the 
fastest experimental rate. As a result, values of T* obtained from simulation must be viewed 
as upper limits to the laboratory glass transition. Velikov, Borick, and AngelF^ estimate 
that simulations overestimate T* by a factor of 1.2-1.6. 

As the temperature of liquid is increased at fixed density, the spinodal temperature 
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Tgpinodai represents the limit of stability of the liquid with respect to formation of a gas. 
Similarly, T% represents the lowest temperature at which a mobile liquid can be consid- 
ered stable with respect to an "ideal" glas d 80 1 81 1 ^. Sastry, Debenedetti, and Stillinger^ 
determined the distributions of "voids", inhomogeneities in the liquid, in simulations of the 
modified Lennard- Jones glass. Voids in the frozen glass appeared only below p* = 0.89, 
which appeared to be the T = extrapolation of the liquid spinodal. Sastry examined this 
phenomena more closely^ and suggested that the spinodal and Kauzmann temperatures 
intersect at T* = 0.16 at p* = 1.08 for a binary Lennard- Jones mixture. 



4. NUMERICAL RESULTS AND COMPARISON TO SIMULATIONS 

T A , T*, T*, and T£ were calculated for several densities p*. We determined the mean 
field stability limit T* A as the temperature at which p A d? = 1.025, where d is the effective 
hard sphere diameter determined from Eqn. [9j The crossover temperature, accounting for 
fractal excitations, T* was determined by the condition that Af = 1.28. T£ was determined 
as the temperature at which the configurational entropy vanishes, i.e., Af = 0. Following 
Lubchenko and Wolynes 48 , the laboratory glass transition temperature T* was calculated 
using Eqn. [5] for the one hour time scale, therefore T* corresponds to Af = 0.82. Figure [T] 
displays T* A as a function of density for the Mulero-Galan-Cuadros equation of state described 
above and for the WCA and KRR separations, as well as for a modification of the KRR 
separation in which A = a(p) at all densities (the reasons for this latter separation are 
discussed below). The results obtained for T A , T*, T* and Tj£ using the other equations of 
state were qualitatively similar to each other. As expected, at the lower densities (p* <~ 
1.01), the KRR and WCA separations give essentially the same results. At higher densities, 
(p* >~ 1-01), however the results from the different potential separations are quantitatively 
different. Because of its discontinuous definition, the KRR separation as it stands would give 
a rather unnatural looking change in slope in the region around (p* ~ 1.01). We therefore 
also explored the predictions of the KRR theory when one chooses A = a(p) at all densities. 
The corresponding results for the KRR separation are shown as solid lines in the figure. All 
three equations of state give T\ ps 1.0 at p* = 1.0. 

It is important to recognize that the stability limit computed from self-consistent phonon 
theory corresponds to a strictly mean field mode coupling limit. Our value of T\ corre- 
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sponds well with that found from detailed mode coupling calculations^. As Reichmann 
and coworkers have pointed outpl, structure based mode coupling calculations give a larger 
value of mode coupling temperature than do fits to simulation data using phenomenological 
mode coupling expressions. The latter phenomenological fits 77 yield a lower temperature 
Tmct = 0.475. This discrepancy between structure based mean field calculations and phe- 
nomenology has been explained by Stevenson, Schmalian and Wolynes^. They have shown 
that fractally shaped excitations allow downhill escape from a mean field minima at temper- 
ature below the mean field T\. Thus, they argue that when fractal or string excitations are 
allowed in random first order transition dynamics, the true stability limit of an aperiodic 
structure is suppressed below its mean field estimate. The crossover to activated dynamics 
in their analysis occurs at a configurational entropy value of S c /N = 1.28, according to 
the simpler estimates based on percolation. We designate this crossover temperature as T* 
and compute it as well as the mean field T\. Its value at p = 0.0 is T* = 0.35. This crossover 
temperature agrees better with the fits of mode coupling phenomenology to simulation. 

Figures |2jj5] display the values of T|, T*, T*, and T^, respectively, predicted by our theory 
using the Mulero-Galan-Cuadros equation of state and both the modified KRR and WCA 
perturbation approaches. The Salzburg- Wood equation of state gives nearly identical values 
for these quantities, while the Carnahan-Starling equation of state gives significantly lower 
values for all the temperatures except for T\. Also shown are values from simulation^ 
for the melting temperature, and the glass transition temperature^ on the simulation 
time scale, T* sim . Examination of our results for T* and T* show that the modified KRR 
separation is in better agreement with simulation than is the WCA separation. Therefore, 
in the remainder of this paper we have used the modified KRR separation. 

The values of T\ are close to the melting temperature at high densities. T* im lies 
between the values of T* and T* as expected. We find T* >sim /T* pa 1.3 close to the estimate 
of AngelP^. The present theory gives T* /T* = 1.5. Stevenson and Wolynes examined^ 
experimental values for T g and confirmed the well-known rule of thumb that T g /T m pa 2/3. 
The present calculations yield lower ratios, T* /T^ pa 0.4 for Zsw and Zmgc and ~ 0.35 for 
Zcs- 

Figure [6] displays T*, T*, and T£ at various values of the pressure. AngelPS has 
tabulated values of Tx/T g from experiment and has found the ratio to typically be ~ 0.78 
for fragile substances at atmospheric pressure. The values from the present theory at p* = 
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are T^fT* = 0.70, 0.72, 0.58, for Zsw, Zmgc, and Zcs, respectively. While the results from 
all equations of state agree at the higher temperatures displayed, the results based on the 
Carnahan-Starling equation of state are substantially different than the other equations of 
state at lower temperatures. The good agreement using equations of state appropriate to 
high density with previous simulations demonstrates the ability of the present theory to 
describe the different temperature regimes in the Lennard- Jones glass. 

The temperature dependences of the predictions for the constant volume and constant 
pressure heat capacity discontinuities ACy and AC P are displayed in Figures [7 10 The 
results found using Z sw (not shown in the figures) and Z M cc are quite similar to each other 
but are distinct (particularly in their high temperature behavior) from the results found 
using Zcs- Table [I] gives the values of ACy/NkB at T* and Table [2] gives the values of 
AC p /NkB at T* . ACp/Nks has been tabulated for a variety of glasses by Stevenson and 
WolynesSSl and the average value of ACp/Nks for glasses of medium fragility at T* per bead 
is 2.85. For the smaller molecular glasses methanol, n-propanol, butyronitrile, ethylene, and 
ethanol, the value per bead is somewhat smaller, averaging 2.1. The p* = values for the 
heat capacity gaps using Zmgc (AC p ~ 1.9) are in agreement with those found for the 
smaller molecular glasses. 

We investigated the conjectures that AC P = k/T and AC P = k! /T 2 that are often used 
in fitting laboratory data. The constants were chosen so that the fits agreed with the 



calculated data at the temperature T*. The results of the fits are shown in Figure 11 It is 
clear that neither of these fitting functions accurately represents the calculated heat capacity, 
though the 1/T* 2 fit is accurate to about 10% above T* except at temperatures close to T\. 
Conversely, the 1/T* fit is better at temperatures below the glass transition temperature 
and can be used determine T£- accurately by integration. The prediction of this procedure 
is shown in the figure. It is clear that the temperature dependence of the calculated heat 
capacities does not precisely follow the simple forms suggested by either gaussian landscapes 
or by elementary two state excitations. 

In discussing the role of intermolecular forces and glasses the empirical correlations that 
have been found between the a-relaxation time r and 1/TV 1 are of keen interest. The 
empirical parameter 7 varies from system to system. The fits to laboratory data give values 
for 7 as small as 0.13 for sorbitol and as large as 8.5 for l,l'-di(4-methoxy-5-methylphenyl)- 
cyclohexane (BMMPC). For a variety of pressures, by varying 7 nearly universal behavior 
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is found when log 10 r is plotted versus l/TV 1 . For a 1/r 12 soft-sphere potential energy, the 
configurational free energy must be a universal function of p* 4 /T by dimensional analysis, 
a fact noticed by several authors^ 2 * 86 * 87 * 88 ! This result is consistent with many of the ex- 
perimental correlations since random first order transition theory provides a link between 
S c and a-relaxation times^. The random first order transition theory expression for the 
relaxation time r = r exp (j^pf^ and = s ^ kg leads to 

! i , ^Nk B 

log 10 r = log 10 t + — (51) 



Figures 12 13 display plots of 1 s 3 versus p* 7 /T using the Mulero, Galan, and Cuadros 
equation of state. The values of 7 shown are best fits to the predicted results. As expected, 
even though the potential is not a power law the behavior is nearly universal. The predicted 
values of 7 do depend on which equation of state is used. For the constant density fits, 
the analysis yields similar values for all 3 equations of state (7 = 5.9 with Zmcg, — 5.8 
with Zsw j an d = 5.9 with Zcs)- There is more variance among the values of 7 when data 
are examined at constant pressure. One finds at constant pressure 7 = 5.8 for the Mulero, 
Galan, and Cuadros equation of state, 7 = 5.5 for the Salzburg- Wood equation of state 
and 7 = 4.7 for the Carnahan-Starling equation of state. The discrepant latter value is 
an indication, seen earlier in our calculations of the characteristic temperatures, that the 
Carnahan-Starling equation of state does not describe the high density glass as well as do 
the other two equations of state. This is not surprising, given the well-known deficiencies 
of this equation of state at high densities. The range of p* 7 /T* displayed in these figures is 
much greater than the range accessible to experiment, due to difficulties in obtaining low 
temperature data. The parameter p*T/T* varies by less than a factor of 2 in experiment; if 
this lmited range were used to display our prediction, the appearance of universality would 
be more striking. The isochoric values of 7 are 10% greater than the corresponding isobaric 
values of 7. 

Casalini and Roland 15 have determined the constant pressure and volume fragilities of 
several glass formers, with the fragility of a glass being defined as 



dlog 10 T 

m = 



(52) 
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Using the random first order transition theory expressions one finds 



UAC v {T g )/Nk B 
(S c (T g )/Nk B y 



m v = — — I „ (53) 



= UZC p{ T s ),N kB 
(S c (T,)/Nk B f 

As always, the temperature T g to be used in these expressions depends on the time scale 
of the experiment. Casalini and Roland evaluate the fragility at a temperature at which 
the relaxation time is 10 s, instead of the more conventional 1 hour definition of T*. We 
have evaluated my and m p on both time scales; Af = 0.82 for the one hour time scale 
and Af = 1.069 for the 10 s time scale. The one hour time scale results are shown in 
Figure 14 In agreement with experiment (see Fig. 3, Huang, et. alP^and Fig. 10, Casalini 



and Roland^, there is little density dependence to m V} while m p is a decreasing function 
of pressure. Direct comparison between theory and laboratory experiment is complicated 
by small errors in the value of D, which we take to be precisely 32, while already the best 
fit to empricial data would differ by about 10%. Also one expects differences between the 
values of AC P found for specific polyatomic substances and pure Lennard-Jonesium. In 
addition, the prediction of absolute fragilities depends quite critically on the corrections 
made to the free energy of the glassy and liquid configurations. We note, however, that 
the ratio m p /m v should be relatively insensitive to all of these effects. Using Z MGC , our 
calculations give m p= o/mv ~ 1.5. On the 10 second time scale, Casalini and Roland fincP^ 
the ratio equal to 1.9 for salol, 1.4 for propylene carbonate, 2.3 for l,l'-di(4-methoxy-5- 
methylphenyl)cyclohexane, 1.9 for phenolphthalein-dimethyl-ether, 1.9 for cresolphthalein- 
dimethyl-ether, and 2.6 for PCB62, a chlorinated biphenyl. 

The Sastry crossing of transitions occurs at low densities where we do not expect the 
Z S w to be accurate. Therefore our search for the Sastry density employed the equations 
of state with Zqs and Zmgc only. We expect the Carnahan-Starling equation of state to 
be the more reliable one at these lower densities. p* spinoda i and T* pinodal were determined by 
locating the temperature at which the pressure as a function of temperature was a minimum 



at a fixed density. Figure [15] displays the spinodal pressure versus density for both equations 
of state and recent molecular dynamics results of Baidakov and ProtsenkcP^. At these low 
densities, Z C s agrees well with simulation, while Z MGC deviates from the simulation data 
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above p* ~ 0.75. The deviation can be traced to the large value of the blip diameter, 
d, that is found at the low spinodal temperatures using Z MGC ; the corresponding values 
of p*d 3 approach the singularities in Zmgc an d render this equation of state unreliable in 



determining the pressures. Figure 16 displays the spinodal and Kauzmann temperatures 
versus density. It is clear that both equations of state predict a crossing of the spinodal 
and Kauzmann lines, at low temperatures, but disagree moderately as to the density of the 
crossing. Z cs predicts a crossing at a density of ~ 0.85 while Z MGC predicts a crossing at 
pa 0.80. Both values are in good agreement with the estimates from simulation®'. 



5. CONCLUSIONS 



Random first order transition theory as instantiated microscopically with self-consistent 
phonon theory and perturbative liquid state theory has been used to predict characteristics of 
the glassy dynamics of a supercoooled Lennard- Jones liquid. The calculations are seen to be 
largely robust with respect to the choice of liquid equation of state and are in good agreement 
with inferences from simulation and experiment. Isochoric and isobaric calculations have 
been performed to distinguish those quantities arising from purely steric interations from 
those due to attractive interactions. We find that the trends of configurational heat capacities 
are reproduced well, but the overall magnitudes are quantitatively 10-20% too small. The 
calculated fragilities follow experimental trends for the density and pressure dependencies 
(demonstrating the present microscopic calculations can accurately reproduce the sign of 
the third derivative of the configurational contributions to the free energy!). Further, the 
spinodal line, as a function of density, is found to intersect the Kauzmann temperature line 
at low densities and temperatures, in agreement with the speculation of Sastry. 

Random first order transition theory is thus able to explain a variety of interrelated 
features of glasses that are directly related to the intermolecular forces, including explaining 
the density-temperature scalings highlighted by Casalini and Roland and the low density 
behavior described by Sastry. The small, but not absent, dependence of the fragility on 
density seen in experiment is reproduced by the theory and the ratio of the isobaric to 
isochoric fragilities is also in reasonably good agreement with experiment. 

The strategy we employ here also provides a route to quantify glassy characteristics of 
complex polyatomic molecular fluids. Molecular liquid theory provides expressions for the 
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equations of state for polyatomic fluids that can be combined with self-consistent phonon 
theory to extend the present treatment to many of the specific molecular liquids commonly 
studied in the laboratory. 
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6. TABLES 



TABLE 1: ACy at T g for the three equations of state used. 





Zsw 


Zmgc 


Zcs 


p* 


a 


AC v /Nk B 


a 


AC v /Nk B 


a 


AC v /Nk B 


1.10 


0.53 


1.47 


0.54 


1.46 


0.38 


0.92 


1.05 


0.42 


1.46 


0.43 


1.45 


0.29 


0.92 


1.00 


0.32 


1.45 


0.33 


1.44 


0.22 


0.90 


0.95 


0.24 


1.43 


0.24 


1.42 


0.16 


0.91 


0.90 


0.17 


1.40 


0.17 


1.40 


0.12 


0.89 


0.85 


0.11 


1.37 


0.12 


1.36 


0.08 


0.88 


0.80 


0.07 


1.30 


0.07 


1.30 


0.05 


0.78 



TABLE 2: AC P at T g using Z SW - 





Zsw 


Zmgc 


Zcs 


p* 


a 


AC p /Nk B 


a 


AC p /Nk B 


1 a 


AC p /Nk B 


25.0 


0.77 


1.76 


0.77 


1.85 


0.64 


1.22 


20.0 


0.68 


1.76 


0.69 


1.86 


0.57 


1.24 


15.0 


0.59 


1.77 


0.59 


1.87 


0.50 


1.25 


10.0 


0.49 


1.78 


0.50 


1.89 


0.42 


1.27 


5.0 


0.39 


1.80 


0.39 


1.93 


0.33 


1.31 


0.0 


0.26 


1.85 


0.27 


2.01 


0.24 


1.41 


-2.0 


0.20 


1.89 


0.21 


2.12 


0.19 


1.52 
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7. FIGURES 
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Figure 2: T\ as a function of density using Zmgc for the modified (A = a(p)) KRR 
(solid line) and WCA (dashed line) separations. T^, the melting temperture^ (filled 
squares) and T* sim , the glass transition temperature at the simulation time scaled 

(filled circles) are also shown. 
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Figure 3: The crossover temperature T* as a function of density using Zmgc for the 
modified (A = a(p)) KRR (solid line) and WCA (dashed line) separations. T^, the 
melting temperture^ (filled squares) and T* sim , the glass transition temperature at the 
simulation time scale^ (filled circles) are also shown. 
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Fi gure 4: T* as a function of density using Z M gc for the modified (A = a(p)) KRR 
(solid line) and WCA (dashed line) separations. T^, the melting temperture^ (filled 
squares) and T* sim , the glass transition temperature at the simulation time scale^ 

(filled circles) are also shown. 
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Figure 5: as a function of density using Zmgc for the modified (A = a(p)) KRR and 
WCA separations. T^, the melting temperture^ (filled squares) and T* sim , the glass 
transition temperature at the simulation time scaled (filled circles) are also shown. 
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Figure 6: T\ (circles), T* (squares), T* (triangles) , and T£ (inverted triangles) versus 
p* predicted using Zsw (solid lines), Zmgc (dashed lines), and Zcs (dash-dot lines). 
The results using Zsw and Zmgc are almost identical. 
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Figure 7: ACy/N ks versus T* using Zmgc at densities p* = 1.10 (circles), p* = 1.00 
(squares), p* = 0.90 (triangles), and p* = 0.80 (inverted triangles). The curves use solid 
(dashed) lines for temperatures below (above) T*. 
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Figure 8: AC V /N k B versus T* using Zcs at densities p* = 1.10 (circles), p* = 1.00 
(squares), p* = 0.90 (triangles), and p* = 0.80 (inverted triangles). The curves use solid 
(dashed) lines for temperatures below (above) T*. 
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Figure 9: AC P versus T* using Zmgc at pressures p* = 25.0 (circles), p* = 15.0 
(squares), p* = 5.0 (triangles), p* = 0.0 (inverted triangles), and p* = —2.0 (diamonds). 
The curves use solid (dashed) lines for temperatures below (above) T*. 
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Figure 10: AC P versus T* using Zqs at pressures p* = 25.0 (circles), p* = 15.0 
(squares), p* = 5.0 (triangles), p* = 0.0 (inverted triangles), and p* = —2.0 (diamonds). 
The curves use solid (dashed) lines for temperatures below (above) T*. 
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Figure 11: Fits of AC P calculated using Zmgc (filled circles) versus 1/T* (dashed line) 
and 1/T* 2 (dash-dot line). The reduced pressure is 0.0 and the fits were performed by 

forcing agreement at the temperature T*. The critical (T* ) and triple (T t * ) 
temperatures from simulation®! are also shown, as are the values of T^, T*, and T£ from 
our calculations. The value of T% obtained by integrating the 1/T* fit to determine the 
temperature at which the configuration entropy would go to zero is also indicated. 
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Figure 12: Plots of l/S'c versus p* 5 - 90 /T* for Zmgc- Values of the density are 1.10 
(circles), 1.05 (squares), 1.00 (triangles), 0.95 (inverted triangles), 0.90 (left-facing 
triangle), 0.85 (right-facing triangle), and 0.80 (diamonds). 
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Figure 13: Plot of 1/S C versus p* 5 - 36 /y* for Zmgc- The pressures are 25.0 (circles), 20.0 
(squares), 15.0 (triangles), 10.0 (inverted triangles), 5.0 (left-facing triangles), 0.0 
(right-facing triangles), and -2.0 (diamonds). 
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Figure 14: Plots my (top panel) and m v (bottom panel) on the one hour time scale 
using Zsw (squares) and Zmgc (triangles). 
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Figure 15: Plot of the spinodal pressures versus density, using Zmgc (circles) and Zcs 
(squares). Also shown are the spinodal pressures from simulationp3 
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Figure 16: Plot of the spinodal (circles, Zmgc, and squares, Zqs) an d Kauzmann 
(inverted triangles, Zmgc, and triangles, Zqs) temperatures versus density. Also shown 

are the spinodal pressures from simulation.^ 
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